Method for vibration analysis

ABSTRACT

Output-only modal testing of an object. Vibrations are excited in said object and measured by a number of vibration sensitive detectors. From the data of the measurements, a spectral density matrix function is determined. From this density matrix, auto spectral densities for the individual modes are identified performing a decomposition based on the Singular Value Decomposition technique. From the auto spectral densities of the individual modes, natural frequencies and damping ratios for the modes can be estimated, and from the singular vectors of the Singular Value Decomposition, the modes shapes can be estimated.

This application claims the benefit of Danish Application No. PA 199901587 filed Nov. 3, 1999 and PCT/DK00/00609 filed Nov. 3, 2000.

BACKGROUND OF THE INVENTION

1. Field of the Invention

The present invention relates to a method for vibration analysisaccording to the preamble of claim 1.

2. Background of the Invention

Modal identification is the process of estimating modal parameters fromvibration measurements obtained from different locations of a structure.The modal parameters of a structure include the mode shapes, natural (orresonance) frequencies and the damping properties of each mode thatinfluence the response of the structure in a frequency range ofinterest.

Modal parameters are important because they describe the inherentdynamic properties of the structure. Since these dynamics properties aredirectly related to the mass and the stiffness, experimentally obtainedmodal parameters provide information about these two physical propertiesof a structure. The modal parameters constitute a unique informationthat can be used for model validation, model updating, quality controland health monitoring.

In traditional modal analysis the modal parameters are found by fittinga model to the Frequency Response Function relating excitation forcesand vibration response. In output-only modal analysis, the modalidentification is performed based on the vibration responses only and adifferent identification strategy has to be used.

Output-only modal testing and analysis is used for civil engineeringstructures and large mechanical structures or structures m operationthat are not easy to excite artificially.

Large civil engineering stares are not easily excited and they are oftenloaded by natural (ambient) loads that cannot easily be controlled ormeasured. Examples of such loads include wave loads on offshorestructures, wind loads on buildings and traffic loads on bridges in suchcases it is an advantage just to measure the natural (or ambient)responses and then estimate the modal parameters by performing anoutput-only modal identification. For civil structures the technique isoften referred to as ambient response testing and ambient responseanalysis.

Application of output-only modal identification instead of traditionalmodal identification gives the user some clearly defined benefits incase of large structures and natural loading. Rather than loading thestructure artificially and considering the natural loading as anunwanted noise source, the latter is used as the loading source. Themain advantages of this technique are:

Testing is less time consuming since equipment for exciting thestructure is not needed.

Testing does not interrupt the operation of the structure.

The measured response is representative of the real operating conditionsof the structure.

When performing output-only modal identification of a structure, theuser can perform the identification in the time domain or in thefrequency domain. For output only identification, the time domaintechniques have been rather dominating since no accurate techniques forfrequency domain identification exists. However, since the frequencydomain supports the physical intuition of the system, i.e. the user canobserve the spectral densities and, thus, directly have an idea of themodes of the system by regarding the spectral peaks, simple and ratherapproximate techniques have been widely accepted for preliminaryanalysis. The most well-known frequency domain technique is theso-called classical approach, also denoted the basic frequency method,or the peak picking method, where the user simply chooses one of thefrequency lines in the spectrum at the appearing peak as the natural(resonance) frequency and then estimates the corresponding mode shape asone of the columns of the spectral density matrix.

The classical approach is based on simple signal processing using theDiscrete Fourier Transform, and is using the fact that well-separatedmodes can be estimated directly from the spectral density matrix at thepeak, as shown in by Julius S. Bendat and Allan G. Piersol in“Engineering Applications of Correlation and Spectral Analysis”, JohnWiley & Sons, 1993.

Other implementations of the technique make use of the coherence betweenchannels as described by A. J. Felber in “Development of a Hybrid BridgeEvaluation System”, Ph.D. thesis, Department of Civil Engineering,University of British Columbia, Vancouver, Canada, 1993. The termchannel is commonly used for the output data of a sensor.

The main advantage of the classical approach compared to otherapproaches, such as two-stage time domain identification technique, orthe one-stage time domain identification techniques, for example theStochastic Subspace Identification algorithms is its user-friendliness.It is fast, simple to use, and gives the user a “feeling” of the data heor she is dealing with.

The two-stage time domain has been described by H. Vold, J. Kundrat, G.T. Rocklin, and R. Russel in “A Multi-Input Modal Estimation AlgorithmFor Mini-Computer”, SAE Technical Paper Series, No. 820194, 1982; by S.R. Ibrahim and E. C. Milkulcik in “The Experimental Determination ofVibration Test Parameters From Time Responses”, The Shock and VibrationBulletin, Vol. 46, No. 5, 1976, pp. 187-196; and by J.-N. Juang and R.S. Papa in “An Eigensystem Realization Algorithm For Modal ParameterIdentification And Modal Reduction”, J. of Guidance, Control andDynamics, Vol. 8, No. 5, 1985, pp. 620-627. The Stochastic SubspaceIdentification algorithms is described by P. Van Overschee and B. DeMoor in “Subspace Identification for Linear Systems”, Kluwer AcademicPublishers, 1996,

The classical technique gives reasonable estimates of naturalfrequencies and mode shapes if the modes are well-separated. However, inthe case of close modes, it can be difficult to detect the close modes,and even in the case where close modes are detected, estimates becomesheavily biased by simple estimation of the mode shapes from one of thecolumns of the spectral matrix.

Further, the frequency estimates are limited by the frequency resolutionof the spectral density estimate, and in all cases, damping estimationis uncertain or impossible.

SUMMARY OF THE INVENTION

It is the purpose of the invention to provide a frequency domain methodfor vibration analysis, i.e. output-only modal analysis, without thesedisadvantages, but where the user-friendliness is preserved.

This purpose is achieved by a method as mentioned by way of introductionand characterised as described in the characterising part of claim 1.

The invention significantly reduces the uncertainty in the estimation ofvibrational modes of an object. Due to its user friendliness and fastobtainable results, the invention is a substantial improvement foroutput-only modal analysis, where the only major difference betweenmodal parameters estimated from traditional modal testing andoutput-only modal analysis is that the output-only modal analysis yieldsunscaled mode shapes.

In the invention, it is assumed that the object has been excited over abroad frequency range by a signal, which has the same intensity at allfrequencies. This kind of excitation is called white noise. As aconsequence of the assumption with regard to the input excitation aswhite noise, the analysis according to the invention is directed to theresponse of the object and uses therefore the well-known term calledoutput-only modal analysis.

A number of techniques are used in the invention. One is the so calledFrequency Domain Decomposition (FDD) which is known from, among others,traditional modal analysis, where the structure is loaded by a knowninput However, in the case of traditional modal analysis the systemmatrix that is decomposed is not the spectral density matrix describingthe responses, but the frequency response function (FRF) matrix relatinginput and output of the system.

As a second technique, the present invention uses the so-called SingularValue Decomposition (SVD) to perform the frequency domain decompositionof the spectral matrix. The SVD is known from numerical mathematics andis used in a number of different applications. However, the mainapplication of this decomposition is to determine the rank of a matrix.In traditional modal analysis, the technique is mainly used to find thenumber of modes, but can also be used to estimate modal properties asdescribed in Shih, C. Y., Y. G. Tsuei, R. J. Allemang and D. L. Brown:“Complex Mode Indicator Function and its Applications to Spatial DomainParameter Estimation”, Mechanical Systems and Signal Proc., Vol. 2, No.4, 1988 and in Shih, C. Y., Y. G. Tsuei, R. J. Allemang and D. L. Brown:“A Frequency Domain Global Parameter Estimation Method for MultipleReference Frequency Response Measurements”, Mechanical Systems andSignal Proc., Vol. 2, No. 4, 1988.

A third technique applied in the analysis of singular value decomposedspectral density functions is the modal assurance criterion (MAC) whichis known from traditional modal analysis and described in Allemang, R.J. and D. L. Brown: “A Correlation Coefficient for Modal VectorAnalysis”, Proc of the 1^(st) International Modal Analysis Conference,IMAC, Orlando, Fla., USA, 1982. In the present invention it is used toisolate auto spectral density functions of the individual modes.

In order to estimate vibrational modes of an object, measurements areperformed in a vibrational excited object using a number of spatiallydistributed vibration sensitive detectors. For example, a building canbe exposed to a broad range of excitation frequencies due to wind ortraffic, justifying the assumption of excitation by white noise.

At suitable predetermined locations, vibration sensitive detectors, forexample accelerometers, are attached to the object and their response ismeasured and transformed to data which are stored in a suitable medium,for example a computer memory.

Those data are subject to further analysis, where the initial step is toachieve a spectrum equivalent to a spectral density function. Thisspectral density function is decomposed into auto spectral densities bya technique involving the modal assurance criterion (MAC), where theauto spectral densities can be interpreted as corresponding toindependent vibrational modes. Then, the auto spectral densities aretransformed from frequency domain to time domain in order to estimatedamping and more accurate natural frequencies.

Though the method itself in practice is relatively simple to use leadingto results in a fast way without the need of cumbersome calculations,the theory behind the method is not obvious and, therefore, one of themajor achievements of the invention. Only a thorough theoreticaltreatment of the mathematical problem behind the theory leads to themethod according to the invention, which is the reason for a detailedexplanation of the theory in the following. After the theoreticalapproach, a practical approach of the method according to the inventionwill be given followed by examples for illustration.

When measuring the response of a structure, the response is sampled, itmeans that the signal is observed at discrete times, each of theobservation times is spaced in time with the sampling interval Δt, andusually a number of observations of the response is obtained in a set ofspatially distributed locations on the structure. When estimating thespectral density of the response, the spectral density is observed in afrequency band from zero frequency to the half of the sample frequencyf_(s)=1/Δt. In this frequency band the response of the structure isinfluenced by the long on the structure and by the dynamic modes in thatfrequency interval. This means, that the spectral matrix of the responseof a structure is a combination of the responses of the modes that ispresent in the observed frequency band. In the following, it is shownthat taking the Singular Value Decomposition (SVD) of the spectralmatrix, the spectral matrix is decomposed into a set of auto spectraldensity functions, each corresponding to a single degree of freedom(SDOF) system. Each SDOF system corresponds to a certain mode of thesystem, i.e., instead of dealing with the more complicated combinationof all modes and all the different spectral densities between theresponses of the different locations of observation on the structure,the problem is reduced to dealing with auto spectral densities of themodal responses of the system. This result is exact in the case wherethe loading is white noise, the structure is lightly damped, and whenthe mode shapes of close modes are geometrically orthogonal. If theseassumptions are not satisfied, the decomposition into SDOF systems isapproximate, but still the results are significantly more accurate thanthe results of the classical frequency domain approach The SDOF autospectral densities are identified using the known modal assurancecriterion (MAC).

Theoretical Background

Obtaining a set of data from one of the detectors attached to the objectof interest yields a set of time related measurements {(y₁,t₁), (y₂,t₂),. . . (y_(n),t_(n))}, where y_(n) is the detector measurement at timet_(n). These measurements are in the following called the systemresponse. The unknown excitation that has generated this system responseis in the following denoted x_(n).

Determining the Spectral Density Function

From this set of measured system response, a frequency dependentspectral density function can be obtained by Fast Fourier Transform(FFT), which is a traditional approach. Assume that the system responsehas been measured at m locations on the object y_(n) is then an m×1vector. Also assume that the Gaussian white noise excitation has beenapplied in r locations, which imply that x_(n) is a r×1 vector.

Due to the assumptions that the unknown excitation x_(n) is Gaussianwhite noise, the m×m spectral density function G_(yy)(iω) of the systemresponse will in the following be defined as

 G _(yy)(iω)=H(iω)G _(xx) H(iω)^(H)−∞<ω<∞  (1)

H(iω) is the m×r Frequency Response Function (FRF) matrix that relatethe applied Gaussian white noise excitation to the system response. Thesuperscript H denotes the Hermitian transpose, i.e. complex conjugateand transpose. The spectral density function G_(xx) of the Gaussianwhite noise excitation is a constant matrix indicating that allfrequencies have been equally excited. The spectral density functionreflects the predominant vibration frequencies ω that occur as peaks inthe spectral density function. Since G_(xx) is a real and symmetricmatrix G_(yy)(iω) will always be Hermitian, i.e. the elements below thediagonal is the complex conjugates of the corresponding elements abovethe diagonal.

However, this method when applied to measured data suffers from theintroduction of the so-called “leakage” biasing in the density spectrumdue to the assumption of periodic data when performing the FFT. Thisbiasing results in a flattening of the before mentioned peaks in thespectral density function. To avoid this biasing, it is possible tocalculate the covariance function for the data set first before a FastFourier Transform. To obtain an unbiased estimate of the covariancefunction however, a tedious calculation is involved.

Therefore, another method for spectral estimation is proposed andapplication of this method together with above mentioned frequencydomain decomposition is a part of the invention. The method that is usedis a combination of traditional FFT and a technique denoted RandomDecrement Transform (RDT). First, a covariance function is estimatedusing the RDT technique, and then, the covariance function istransformed to frequency domain by the FFT. Since the covariancefunction is a decaying function in time, it is well known that byperforming FFT on this covariance function instead of the data itself,the leakage bias of the resulting spectral density function will besignificantly decreased. In fact, if the covariance function hasdecreased to zero at half of the maximum time lag, then the leakage biasis exactly zero.

The advantages of using the RDT technique are:

Instead of performing time consuming multiplication, as is necessary fortraditional covariance function estimation, the RDT uses only additionswhich are much faster calculations for normal computers;

Instead of being constrained by being forced to use all the signal inone way only, the user can weight the different parts of the time signalby using different kinds of trigger level, where the trigger level isexplained in the following below;

Instead of getting only one estimate of the covariance function, theuser gets two independent estimates, one for the covariance functionitself, and one for its derivative.

It has to be mentioned here that modes at higher frequencies areweighted more, if the covariance function in the above method isreplaced by the derivative of the covariance function. Therefore, insome cases, using the derivative of the covariance function can be anadvantage.

When the user of the RDT technique is to perform an estimate of thecovariance function, he has to choose a so-called trigger criterion. Thetrigger criterion can for instance be that the signal is within acertain intensity range. Then, if the signal is within this range, thetechnique selects the data around the time where the trigger criterionis satisfied, and incorporates these data in the estimation process.Thus, the user has the possibility to estimate one covariance functionfor low amplitudes and one for large amplitudes. This is of importancein cases where the user might want to investigate whether the structureunder consideration has been subject to vibrations with such largeamplitudes, that the structural response has been become non-linear,i.e. parts of the structure has been yielding (steel) or cracking(concrete or masonry). In this case, if the structural response isnon-linear, the user will observe a difference of the spectral densityfunction estimates for the low and for the high amplitudes.

The RDT technique is a well known technique for estimation of systemtime functions. However, it is not common to transform the timefunctions to frequency domain, and the above mentioned advantages bydoing so are not known.

Once, the spectral density function has been obtained the followingtools are applied.

Theoretical Background of Frequency Domain Decomposition

Assume that the spectral density function G_(yy)(iω) of the systemresponse has been estimated at discrete frequencies ω=ω_(j), forω_(j)≧0. G_(yy)(iω) is then described by a set of spectral densitymatrices, one for each of the discrete frequencies.

Consider the Hermitian spectral density matrix for ω=ω_(j). This matrixis then decomposed by taking the Singular Value Decomposition (SVD) ofthe matrix

G _(yy)(iω _(j))=U _(j) S _(j) U _(j) ^(H) =s _(j1) u _(j1) u _(j1) ^(H)+s _(j2) u _(j2) u _(j2) ^(H) +, . . . , +s _(jm) u _(jm) u _(jm)^(H)0≦ω≦∞  (2)

where the matrix U_(j)=└u_(j1),u_(j2), . . . , u_(jm)┘ is a unitarymatrix holding the orthogonal singular vectors u_(jk) as its columns,and S_(j) is a diagonal matrix holding the scalar singular valuess_(jk). The singular values and corresponding vectors are always sortedso that the s_(j1) always will be the largest. This decomposition isperformed for all the estimated spectral density matrices.

Near a peak, corresponding to one mode, the Eigenvalue of that modebecomes dominating resulting in a loss of rank of the spectral matrix atthat particular frequency. This loss of rank is detected by the SVD by asimilar loss of singular values having values significantly differentfrom zero. If only one mode is completely dominating at a specificfrequency ω=ω_(j), only s_(j1) will be significantly different fromzero. As shown below, an estimate of the mode shape of this dominatingmode is given by the corresponding singular vector u_(j1) assuming thatthe mode is lightly damped. The mode shape is the observable part of theEigenvector that is associated with the dominating Eigenvalue.

If two close modes are dominating at ω=ω_(j) then s_(j1) and s_(j2) willboth be significant different from zero and corresponding mode shapesare given by the singular vectors u_(j1) and u_(j2) assuming that themode shapes are orthogonal as u_(j1) and u_(j2) by definition are. Thisapproach can in principle be extended to as many close modes as thereare measurement channels m.

It is only at the peak of a dominating mode that u_(j1) will be the modeshape estimate. When moving away from the peak, the influence of theEigenvalues of the other modes start to interfere. This means that thesingular vector that is an estimate of the mode shape now can be any oneof the singular vectors. However, by using the so-called Modal AssuranceCriterion (MAC) it is possible to track which of the singular vectorsthat currently matches the u_(j1) at the peak. From this tracking it isthen possible to determine which singular values that relates to themode at the peak, and in this way construct the auto spectral densityfunction of the particular mode. Due to the presence of noise it mightnot be possible to track the corresponding singular vectors over thehole frequency range. In this case the untracked part of the autospectral density function is simply set to zero.

That statement that the first singular vector is an estimate for themode shape can be verified by the following proof

Let “−” and superscript T denote complex conjugate and transpose,respectively. The FRF can then be written in partial fractions, i.e.pole/residue form $\begin{matrix}{{H\left( {\quad \omega} \right)} = {{{\sum\limits_{k = 1}^{n}\frac{R_{k}}{{\quad \omega} - \lambda_{k}}} + \frac{{\overset{\_}{R}}_{k}}{{\quad \omega} - {\overset{\_}{\lambda}}_{k}}\quad - \infty} \leq \omega \leq \infty}} & (3)\end{matrix}$

where n is the number of modes. λ_(k) is the Eigenvalue (pole) and R_(k)is the m×r residue matrix

 R _(k)=φ_(k)γ_(k) ^(T)  (4)

where φ_(k),γ_(k) is the mode shape vector and the modal participationvector, respectively. By inserting (3) into (1) the following relationis obtained $\begin{matrix}{{G_{yy}\left( {\quad \omega} \right)} = {{{\sum\limits_{k = 1}^{n}{\sum\limits_{s = 1}^{n}{\left\lbrack {\frac{R_{k}}{{\quad \omega} - \lambda_{k}} + \frac{{\overset{\_}{R}}_{k}}{{\quad \omega} - {\overset{\_}{\lambda}}_{k}}} \right\rbrack {G_{xx}\left\lbrack {\frac{R_{s}}{{\quad \omega} - \lambda_{s}} + \frac{{\overset{\_}{R}}_{s}}{{\quad \omega} - {\overset{\_}{\lambda}}_{s}}} \right\rbrack}^{H}}}}\quad - \infty} \leq \omega \leq \infty}} & (5)\end{matrix}$

By multiplying the two partial fraction factors, and making use of theHeaviside partial fraction theorem, the spectral density functionG_(yy)(iω) can, after some mathematical manipulations, be reduced to apole/residue form as follows $\begin{matrix}{{G_{yy}\left( {\quad \omega} \right)} = {{{\sum\limits_{k = 1}^{n}\frac{A_{k}}{{\quad \omega} - \lambda_{k}}} + \frac{{\overset{\_}{A}}_{k}}{{\quad \omega} - {\overset{\_}{\lambda}}_{k}} + \frac{B_{k}}{{- {\omega}} - \lambda_{k}} + \frac{{\overset{\_}{B}}_{k}}{{{- }\quad \omega} - {\overset{\_}{\lambda}}_{k}}\quad - \infty} \leq \omega \leq \infty}} & (6)\end{matrix}$

where A_(k) and B_(k) are the k'th residue matrices of the spectraldensity function G_(yy)(iω). The contribution from the B_(k) elements ismuch smaller than the contribution from the A_(k) elements and will inthe following be regarded as negligible.

As the spectral density function G_(yy)(iω) itself the residue matricesare m×m Hermitian matrices.

The residue matrix A_(k) is given by $\begin{matrix}{A_{k} = {R_{k}{G_{xx}\left( {{\sum\limits_{s = 1}^{n}\frac{R_{s}}{{- {\overset{\_}{\lambda}}_{k}} - \lambda_{s}}} + \frac{{\overset{\_}{R}}_{s}}{{- {\overset{\_}{\lambda}}_{k}} - {\overset{\_}{\lambda}}_{s}}} \right)}^{H}}} & (7)\end{matrix}$

The contribution to the residue from the k'th mode is given by$\begin{matrix}{A_{k} = \frac{R_{k}G_{xx}R_{k}^{H}}{2\left( {2\quad \pi \quad f_{k}\zeta_{k}} \right)}} & (8)\end{matrix}$

where the denominator of (8) is two times minus the teal part of theEigenvalue λ_(k)=−2πf_(k)ζ_(k)+i2πf_(k) {square root over (1−ζ_(k) ²)}with f _(k),ζ_(k) being the natural frequency and damping ratio,respectively. As it appears, this term becomes dominating when thedamping is light, i.e. when ζ_(k) tends to zero. Thus, is case of fightdamping, the residue becomes proportional to the mode shape vector

A _(k) ∝R _(k) G _(xx) R _(k) ^(H)=φ_(k)γ_(k) ^(T) G _(xx){overscore(γ)}_(k)φ_(k) ^(H) =d _(k)φ_(k)φ_(k) ^(H)  (9)

where d_(k) is a scalar constant that is non-negative since G_(xx)always will be at least semi-definite.

At a certain frequency ω, only a limited number of modes will contributesignificantly, typically one or two modes. Let this set of modes bedenoted by Sub(ω). Thus, in the case of a lightly damped object, thespectral density function G_(yy)(iω) can always be written$\begin{matrix}{{G_{yy}\left( {\quad \omega} \right)} = {{{\sum\limits_{k \in {{Sub}{(\omega)}}}\frac{d_{k}\phi_{k}\phi_{k}^{H}}{{\quad \omega} - \lambda_{k}}} + \frac{{\overset{\_}{d}}_{k}{\overset{\_}{\phi}}_{k}\phi_{k}^{T}}{{\quad \omega} - {\overset{\_}{\lambda}}_{k}}\quad - \infty} \leq \omega \leq \infty}} & (10)\end{matrix}$

Now (10) describes the spectral density function from −∞ to ∞. However,in practice, only the positive part of the spectral density function isconsidered, which reduces (10) to $\begin{matrix}{{G_{yy}\left( {\quad \omega} \right)} = {{\sum\limits_{k \in {{Sub}{(\omega)}}}\frac{d_{k}\phi_{k}\phi_{k}^{H}}{{\quad \omega} - \lambda_{k}}} = {{\sum\limits_{k \in {{Sub}{(\omega)}}}{\frac{d_{k}}{{\quad \omega} - \lambda_{k}}\phi_{k}\phi_{k}^{H}\quad 0}} \leq \omega \leq \infty}}} & (11)\end{matrix}$

This is a modal decomposition of the spectral matrix. The expression issimilar to the results one would get directly from (2) under theassumption of independent white noise input, i.e. a diagonal spectralinput matrix. For each frequency and each eigenvalue, value, the ratioin front of the mode shape product in (11) will be a positive constantthat can always be made real by proper scaling of the corresponding modeshapes.

Identification Algorithm

Near a peak corresponding to the k th mode in the spectrum say atω=ω_(j), this mode or maybe a possible close mode will be dominating. Ifonly the k th mode is dominating, there will only be one term in (11).Thus, in this case, the first singular vector u_(j1) is an estimate ofthe mode shape

{circumflex over (φ)}=u _(j1)  (12)

and the corresponding singular value is the current value of auto powerspectral density function of the corresponding single degree of freedomsystem at that specific frequency, refer to (11).

Define the MAC value between two vectors as $\begin{matrix}{{{MAC}\left( {\hat{\phi},u_{nm}} \right)} = \frac{{{\hat{\phi}}^{H}u_{nm}}}{\sqrt{{\hat{\phi}}^{H}\hat{\phi}}\sqrt{u_{nm}^{H}u_{nm}}}} & (13)\end{matrix}$

This value describes the correlation between the two vectors. If thevectors are similar except for a constant scaling the MAC value is 1. Ifthey are complete orthogonal the value will be 0.

Using the MAC value in (13) the rest of the auto spectral densityfunction is identified around the peak by comparing the mode shapeestimate {circumflex over (φ)} in (12) with the singular vectors for thefrequency lines around the peak. As long as a singular vector u_(nm) isfound that has MAC value near 1 with {circumflex over (φ)}, thecorresponding singular value belongs to the auto spectral densityfunction of that specific mode

If at a certain distance firm the peak of the mode none of the singularvalues has a singular vector with a MAC value larger than a certainlimit value Ω, the search for matching parts of the auto spectraldensity function is terminated. The remaining unidentified part of theauto spectral density function is set to zero. From this fully orpartially identified auto spectral density function of the SDOF system,the natural frequency and the damping are obtained by taking the autospectral density function back to time domain by inverse discreteFourier transform.

From the free decay time domain function, which is also the autocorrelation function of the SDOF system of the k th mode, the naturalfrequency and the damping ratio is found by estimating crossing timesand logarithmic decrement First, all extremes r_(n), both peaks andvalleys, on the auto correlation function are found. The logarithmicdecrement δ is then given by $\begin{matrix}{\delta = {\frac{2}{n}\ln \quad \left( \frac{r_{0}}{r_{n}} \right)}} & (14)\end{matrix}$

where r₀ is the initial value of the auto correlation function and r_(n)is the n th extreme. Thus, the logarithmic decrement and the initialvalue of the correlation function can be found by linear regression onnδ and 2ln(|r_(n)|), and the damping ratio is given by the well knownformula $\begin{matrix}{\zeta_{k} = \frac{\delta}{\sqrt{\delta^{2} + {4\quad \pi^{2}}}}} & (15)\end{matrix}$

A similar procedure is adopted for determination of the naturalfrequency. The natural frequency is found by making a linear regressionon the crossing times and the times corresponding to the extremes andusing that the damped natural frequency f_(k) ^(d) and the undampednatural frequency f_(k) is related by $\begin{matrix}{f_{k} = \frac{f_{k}^{d}}{\sqrt{1 - \zeta_{k}^{2}}}} & (16)\end{matrix}$

The extreme values and the corresponding times are found by quadraticinterpolation, whereas the crossing times where found by linearinterpolation.

BRIEF DESCRIPTION OF THE DRAWINGS

Having described the theory behind the method according to theinvention, the invention will be explained further in detail withreference to the following drawing, where

FIG. 1 shows a density spectrum with two close modes,

FIG. 2 shows a partial identification of the two peaks from FIG. 1,

FIG. 3 shows a linear regression on extremes for estimation of damping(top), the time domain free decay obtained by inverse FFT, and anestimated damping envelope (bottom),

FIG. 4 illustrates case 2 with closely spaced modes. Top: Partialidentification of SDOF auto spectral density. Bottom: Corresponding freedecay with damping envelope.

FIG. 5 illustrates case 3 with closely spaced modes, but where only avery limited part of the SDOF density is identified,

FIG. 6 illustrates case 4 with moderately spaced modes; mode shapes notorthogonal,

FIG. 7 illustrates case 5 with moderately spaced modes; correlatedinput,

FIG. 8 illustrates the different vibrational modes for a building.

DETAILED DESCRIPTION OF THE PREFERRED EMBODIMENTS

The technique according to the invention and as described above isillustrated on a case with two closely spaced modes in the densityspectrum. The response of a two-degree-of-freedom system is simulatedusing a vector ARMA model as described by P. Andersen, R. Brincker, andP. H. Kirkegaard in “Theory of Covariance Equivalent ARMAV Models ofCivil Engineering Structures”, Proc of the 14^(th) International ModalAnalysis Conference, IMAC, Dearborn, 1996. Furthermore, it is assumedthat these two degrees of freedom are loaded by Gaussian distributedwhite noise uncorrelated processes. The output of the simulation is arealisation of system response in two channels, i.e. theoretical datasets from only two detectors are regarded.

The first case considered is a case with a reasonable spacing betweenthe two modes. The spectral density function 2 of the system response ofthe first of the two channels is shown in the top spectrum of FIG. 1,and the singular values 4, 6 of the decomposed spectral densityfunctions, which corresponds to the spectral density matrix, dependenton frequency are shown in the bottom spectrum of FIG. 1. As it appears,the two modes 8, 10 are clearly visible as peaks in both spectra Partialidentifying of the auto spectra) density functions of the two SDOFsystems using the MAC as described above yields the result as shown inFIG. 2, where each of the peaks 8, 10 are identified. The two singlemodes of the system can now be identified by the auto spectral densityfunctions 12 and 14. These auto spectral density functions 12 and 14 asillustrated in FIG. 1 approximate to a very high degree the moretheoretical SDOF auto spectral density functions and will in thefollowing be treated as being the SDOF auto spectral density functions.The first single mode illustrated by the first auto spectral densityfunction 12 corresponds approximately to the combination of that part12′ of the first singular value 4 which is to the left the line 16 at 15Hz and that part 12″ of the second singular value 6 which is to theright of the line 16 at 15 Hz. In analogy, the second single modeillustrated by the second auto spectral density function 14 correspondsapproximately to the right part 14′ of the first singular value 4 andthe left part 14″ of the second singular value 6.

Taking the inverse discrete Fourier transform of the above partiallyidentified SDOF auto spectral density function 12 of the first modeyields the corresponding auto correlation function estimate 18 as shownin FIG. 3, bottom. Top part of the same figure shows the linearregression on the extremes. The straight line 20 is the best fit line ina least squares sense. The so-called damping envelope 22 as shown in thebottom part of the figure is obtained by an exponential calculation withthe values of the straight line 20 as arguments. The estimated dampingenvelope 22 is congruent with the positive extreme values of the autocorrelation function 18.

As it appears, the procedure is quite straight forward and the user hasa clear impression of the validity of the estimation simply byinspecting the plots.

Exact and identified natural frequencies and damping ratios are shown inTables 1 and 2 for this first case and for the four other cases asdescribed in the following. f₁ and f₂ are the frequencies for the firstand the second single modes with corresponding damping ratios ζ₁ and ζ₂.

TABLE 1 Exact and estimated natural frequencies Exact Exact EstimateEstimate Case ƒ₁ (Hz) ƒ₂ (Hz) {circumflex over (ƒ)}₁ (Hz) {circumflexover (ƒ)}₂ (Hz) 1 14.235 15.916 14.241 15.907 2 15.532 15.916 15.50815.896 3 15.532 15.916 15.526 15.888 4 14.309 15.834 14.433 15.679 514.235 15.916 14.226 15.931

TABLE 2 Exact and estimated damping ratios. Exact Exact EstimateEstimate Case ζ₁ (%) ζ₂ (%) {circumflex over (ζ)}₁ (%) {circumflex over(ζ)}₂ (%) 1 0.894 1.000 0.913 1.164 2 0.976 1.000 0.966 0.938 3 0.9761.000 0.670 0.579 4 2.225 2.539 2.589 3.113 5 0.894 1.000 1.083 0.986

The second case considered is the case of closely spaced modes as shownin FIG. 4. In this case, it is assumed that a reasonable part of theSDOF auto spectral density function can be identified on both sides ofthe peak of the considered mode. This is possible in the most cases byspecifying a lower Ω-value as the MAC limit. In this case, theidentification is also straight forward and the identified values offrequency and damping compares well with the exact values, see Tables 1and 2.

In the third case, it is assumed that only a quite small part of theSDOF auto spectral density function can be estimated. This can be thecase if the spectral density is noisy due to limited data, or if noiseis contaminating the signal. In this case, however, since the data aresimulated data meeting all basic assumptions of the technique, theidentified SDOF density function shown in FIG. 5 was obtained by using arather high value of the MAC limitΩ. As it appears, since the number ofactive frequency lines in the spectrum is cut significantly, the freedecay in the time domain becomes truncated to a degree where the dampingbecomes underestimated, Table 2.

The fourth case, as shown in FIG. 6, illustrates the influence ofnon-orthogonal modes. In theory, to give exact results, the presentmethod requires that the modes are orthogonal. All other casesconsidered in this paper have orthogonal modes. For the modes consideredin case four, the MAC matrix is ${MAC} = \begin{bmatrix}1.0000 & 0.4226 \\0.4226 & 1.0000\end{bmatrix}$

In this case, the SVD still splits the spectral matrix into orthogonalcomponents. This means, that even though the dominant singular value andthe corresponding singular vector is a good estimate of the modalproperties, the second singular value and the corresponding vector arenot so closely related to the physics of the system. Thus, the rightmost part 24 of the SDOF auto spectral density function of the left mode26 is badly estimated. Even though this is the case, the estimates offrequency and damping are still close to the exact values, Tables 1 and2.

For the last considered case, case five, the loading is moderatelycorrelated. In case of correlated input the modal decompositionperformed by the method is approximate. In most practical cases however,like wind loads, wave loads or traffic loads, it is known that amoderate spatial correlation is present. Thus it is important to knowthe amount of influence such correlation might have on the results. Inthis case the correlation matrix between the two stochastic processesloading the system was $C = \begin{bmatrix}1.0000 & 0.4724 \\0.4724 & 1.0000\end{bmatrix}$

The results of the estimation of the damping ratio of the first mode areshown in FIG. 7. Again, a certain distortion of the identified autospectral density of the associated SDOF system in the overlapping regionbetween the two modal peaks is seen. However, the influence is rathersmall, the estimation of frequency and damping is straight forward, andthe estimated values are close to the exact values, Tables 1 and 2.Thus, moderate correlation does not seem to significantly influence thequality of the results.

For illustration, the first five vibrational modes are shown in FIG. 8,where FIG. 8a shows the building 80 without vibrations, FIG. 8b showsthe same building 80 with an overlay 81 illustrating the first vibrationmode, FIG. 8c shows the same building 80 with an overlay 82 illustratingthe second vibration mode, FIG. 8d shows the same building 80 with anoverlay 83 illustrating the third vibration mode, FIG. 8e shows the samebuilding 80 with an overlay 84 illustrating the fourth vibration mode,and FIG. 8f shows the same building 80 with an overlay 85 illustratingthe fifth vibration mode.

Application of the Invention in Civil and Mechanical Engineering

The typical reason for applying the invention is to obtain informationabout the stiffness or mass of a structure.

The modal parameters are closely related to the mass and to thestiffness of a structure. The stiffer the structure, the higher thenatural frequencies, and the heavier the structure, the lower thenatural frequencies. The mode shapes reflect the mass and stiffnessdistribution, and the damping ratio reflects the internal dampingproperties of the structure.

Stiffness information can also be obtained for instance by measuring therelation between deformation and static loading, but the advantage ofusing the dynamic response is that the testing is much easier and inmany cases the stiffness information obtained in this way is moreaccurate.

Once the stiffness of the structure has been obtained by the invention,one of the most important applications of this kind of stiffnessinformation is for validation of numerical models and for qualitycontrol of parts produced in large numbers.

It is common to model the dynamic behaviour of structures by usingnumerical models like finite Element models. This is common practice inmechanical as well as in civil engineering. However, no matter how welldeveloped the modelling tools are, still some uncertainties existconcerning material behaviour and structural joints. For this reason, itis of great practical value to validate these models by comparing themodelled modal parameters with the modal parameters obtained fromexperiments. If the model compares well, then the model is consideredreliable, and can provide a good basis for simulating the behaviour ofthe structure under different loading conditions.

Examples of this kind of applications are typically expensive prototypeslike satellites, rockets, aeroplanes, and large bridges, buildings anddams.

If a structural part is produces in large numbers, then it is often ofgreat value to have a tool to check if every individual part has thesame properties as the rest of the production. In this case the modalproperties of the individual part are compared to the average modalproperties of the production, and if they compare well, then the qualityof the individual part has been documented. This kind of quality controlis used on all kinds of production where the stiffness properties areimportant, or where minor flaws can be catastrophic, examples are woodmembers, concrete members, brake disks, engine blocks etc.

Stiffness information obtained by the invention can also be used formodel updating. If for some reason a numerical model does not comparewell with the measured modal parameters, it might be useful to improvethe model so that the modal parameters of the model are closer to theobserved parameters. This is especially the case if the model is to beused for response simulation.

This kind of model improvement is called model updating. Typicalexamples of this kind of model improvement is modelling of large enginesand turbines, modelling of aeroplanes, rockets and satellites, andmodelling of large civil engineering structures like large bridges,buildings and dams.

Until recently this kind of model improvement was only possible usingtraditional modal testing, i.e., when the structure can be loaded byforces that can be measured and controlled. However, since output-onlymodal analysis now is available, the output only modal identificationcan be performed on the large civil engineering structures since themodal parameters can be estimated from the measured responses only, andthus, better models can be achieved for these structures. This meansthat the response of these structures to dangerous natural loads such asstorms and earth quakes can be modelled better and money can be saved orlarger safety can be documented.

Stiffness information obtained by the invention might also be used forhealth monitoring of structures. Health monitoring is of importancewhenever an expensive structure has to operate for a long time and whereits mechanical health can be uncertain due to fatigue, wear orenvironmental influence such as frost and chemical reactions.

In such cases, it is of large importance to know when the health of thestructure has decreased, and also to know when the health is so muchdecreased that the right moment for repair has appeared.

Nowadays, most of these health monitoring tasks are performed by visualor manual inspection. However, since the measuring equipment is becomingcheaper and more accurate, and since the facilities for transmittingthis kind of information is already developed it is anticipated that themain part of this kind of monitoring can be automated.

A lot of research is going on in this area, and these years, severalinternational conferences concern this technology. In the researchcommunity, it is commonly accepted, that one of the corner stones insuch systems will be the automated modal identification using outputonly techniques.

Examples of applications can be large engines and turbines, aeroplanes,rockets, satellites and all kinds of large expensive structures likelarge bridges, buildings and dams. However in the long run also smallerstructures like highways bridges might be monitored using this kind oftechnology. This opens a very large market for output onlyidentification.

An important application could be the evaluation of structural safetyafter an earth quake where possible damages could exist. Then making ahealth evaluation immediately after the earth quake it can be documentedthat no significant structural changes has occurred, and the buildingcan be authorised for use right away.

In case it has been detected, that the modal parameters have changedsignificantly, then it could be analysed what the reasons are for thesechanges. This is called damage detection.

If a numerical model has been calibrated (updated as described above)for the structure in its virgin state (when the structure was new), thismodel can be used for making systematic changes of the numerical modelto see where in the structure possible stiffness changes can explain theobserved changes of the modal parameters.

If a clear indication of the location of the damage can be made,substantial effort in finding the damage can be saved. This isespecially true for structures where the main part of the structure ishidden like it is the case for office and apartment buildings.

If not only the location can be estimated, but also the type and thesize of the damage can be estimated by performing these numericalsimulations, the consequences for the safety of the structure can beestimated. This is of great importance for large civil engineeringstructures in earth quake active areas and areas with typhoons andHurricanes. Thus, if this becomes possible, then immediately after anearth quake, the safety of the structure can be assessed, and the costof the repair can be estimated.

Application of the Invention in Other Fields

The invention can be successfully used for analyses of any data that isthe response of a system where the loading is unknown.

Any mechanical system can be analysed, it can be typical mechanical orcivil engineering structures as described above, but can also begeological sites or parts of rock or earth, bio-mechanical systems, likehuman bodies, animal bodies and plants and parts hereof.

It can be used for analysis of the radiation of information that dependson the movements of any mechanical systems for instance radiation fromplanets and stars.

The invention can also be used for the analyses the response ofmicro-particles like atoms and molecules.

There are many fields in electronic engineering where the techniquecould be useful. In this engineering field identification has beenextensively used, but emphasis has been on techniques where the inputwas known, on time-domain techniques, and on cases where the response isobserved by only one signal (single output). In case of multiple outputobservations, the invention could have many applications in theelectronic engineering fields in control and identification of systems.In particular, the invention could be used for identification ofacoustic modes, i.e. the resonance frequencies of a room.

In economics the fluctuations of stocks and interest rates are oftendescribed by stochastic dynamic models. However, like for electronicengineering, the techniques used in economics are usually time domaintechniques, and usually limited to one-output cases. For analysis oflarge sets of economic data, and especially for cases where severalobservation have been made (multi output), the invention could be apowerful tool for analysis of such data.

Implementation of the Invention

The method according to the invention is well suited to be implementedin a vibration analysing computer program, where the measured data fromthe object serve as an input and where the vibrational modal functionsare calculated together with damping factors related to the single modesand corresponding natural frequencies.

What is claimed is:
 1. Method for testing an object with respect to itsvibration modes, comprising exciting vibrations in said object,measuring said excited vibrations by at least one vibration sensitivedetector, and obtaining data quantitative of said excited vibrations ofsaid object, wherein the method further comprises; determining aspectral density function identifying in said spectral density functionat least one auto spectral density function corresponding to a singledegree of freedom system and obtaining a decomposition of said excitedvibrations into mutually independent vibration modes.
 2. Methodaccording to claim 1, wherein the decomposition comprises the SingularValue Decomposition of the spectral density matrix at individualfrequency lines.
 3. Method according to claim 1, wherein said autospectral densities are identified using the known modal assurancecriterion (MAC).
 4. Method according to claim 1, wherein said autospectral function is transformed by Inverse Fourier Transform fordetermination of the vibration damping characteristics of said objectand/or the natural frequency of the vibration mode.
 5. Method accordingto claim 1, wherein the determination of said spectral density functioncomprises Fast Fourier Transform.
 6. Computer program using a methodaccording to claim
 1. 7. A method according to claim 1, wherein saidobject is biomechanical systems as the human body, animal body, plants,or parts thereof, civil engineering structures such as buildings,bridges or dams, structural parts such as wood members, concretemembers, or metal members, or members in mechanical structures likeheavy machinery, rotating machinery, brake disks, engine blocks,turbines, aeroplanes, rockets, or satellites.
 8. Method for testing anobject with respect to its vibration modes, comprising excitingvibrations in said object, measuring said excited vibrations by at leastone vibration sensitive detector, and obtaining data quantitative ofsaid excited vibrations of said object, wherein the method furthercomprises determining a spectral density function identifying in saidspectral density function at least one auto spectral density functioncorresponding to a single degree of freedom system and obtaining adecomposition of said excited vibrations into mutually independentvibration modes, wherein the determination of said spectral densityfunction comprises Fast Fourier Transform, and wherein the determinationof said spectral density function further comprises Random DecrementTransform of said data preceded by the Fast Fourier Transform.
 9. Methodaccording to the claim 8, wherein said the determination of saidspectral density function comprises Random Decrement Transform atdifferent trigger levels as a control for non-linearities of the systemresponse.
 10. Method for testing an object with respect to its vibrationmodes, comprising exciting vibrations in said object, measuring saidexcited vibrations by at least one vibration sensitive detector, andobtaining data quantitative of said excited vibrations of said object,wherein the method further comprises determining a spectral densityfunction identifying in said spectral density function at least one autospectral density function corresponding to a single degree of freedomsystem and obtaining a decomposition of said excited vibrations intomutually independent vibration modes, wherein said determination of saidspectral density function comprises the derivative of the covariancefunction obtained by the Random Decrement Transform of the time data.